##################################################################
# Nadiya Kostyuk
# Running Analysis
# R version: 4.1
# 11/19/25
##################################################################

rm(list=ls())

## Install & load packages (all at once)
list.of.packages <- c('dplyr', 'MASS', 'lmtest', 'stargazer', 'stringr', 
                      'lubridate', 'ggplot2', 
                      "tidyr","foreign","survival","survminer","broom","knitr",
                      "corrplot",'tidyverse','forcats',
                      'reshape2', 'lme4')
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]; if(length(new.packages)){install.packages(new.packages,dependencies=TRUE)}
lapply(list.of.packages, require, character.only = TRUE)


# setting directory: 
setwd("")


# loading data
load('data/data_final.RData')


########################
# Table 1:
########################


##########################
# Running Analysis 
# For Contested Territories
# Maneuver versus Attrition
# H1A-H1C
###########################

# subsetting data to contested territories:
data_d_contested_80 = data_final %>% 
  filter(territorial_control_80_1dlag == 'Contested')
dim(data_d_contested_80)


mod_list = list()

# basic model: 
mod1a = glmer(
  form = disruption_intentional_3indicators  ~
    disruption_intentional_3indicators_1dlag
  + as.factor(war_phases_1dlag)
  + (1|ADM1_NAME), 
  data = data_d_contested_80, 
  family = binomial()
)
summary(mod1a)
mod_list[[1]] = mod1a

# model with control
mod1b = glmer(
  form = disruption_intentional_3indicators  ~
    disruption_intentional_3indicators_1dlag
  + as.factor(war_phases_1dlag)
  + log_rus_mil_t_1dlag
  + holiday_dummy
  + as.factor(WEEKEND)
  + (1|ADM1_NAME), 
  data = data_d_contested_80, 
  family = binomial()
)
summary(mod1b)
mod_list[[2]] = mod1b


##########################
# Running Analysis 
# For Contested versus Ukrainian-controlled Territories
# Punishment versus Maneuver/Attrition
# H2
###########################

# excluding Russia-occupied territories:
data_d_no_ruscontrol_80 = data_final %>% 
  filter(territorial_control_80_1dlag != 'Russian occupation')
dim(data_d_no_ruscontrol_80)

# basic model: 
mod1c = glmer(
  form = disruption_intentional_3indicators  ~
    disruption_intentional_3indicators_1dlag
  + as.factor(territorial_control_80_1dlag)
  + (1|ADM1_NAME), 
  data = data_d_no_ruscontrol_80, 
  family = binomial()
)
summary(mod1c)
mod_list[[3]] = mod1c


# model with control
mod1d = glmer(
  form = disruption_intentional_3indicators  ~
    disruption_intentional_3indicators_1dlag
  + as.factor(territorial_control_80_1dlag)
  + log_rus_mil_t_1dlag
  + holiday_dummy
  + as.factor(WEEKEND)
  + (1|ADM1_NAME), 
  data = data_d_no_ruscontrol_80, 
  family = binomial()
)
summary(mod1d)
mod_list[[4]] = mod1d


# saving model results: 
#saveRDS(mod_list, file = 'analysis/output/tables/main_rus.RDS')

################################
# Code for generate Table 1 for
# the Online Appendix:
stargazer(mod1a, mod1b, mod1c, mod1d)


############################################
# Turning results from above into Table 1 in the main manuscript: 

# # loading functions:
source('analysis/code/step3_create_latex_table_function.R')

# loading results:
rez = readRDS('analysis/output/tables/main_rus.RDS')
mod_list = rez
mod_list

mod_res_wide = create_latex_table_glm(mod_list)
mod_res_wide

# saving the model:
#writeLines(capture.output(xtable(mod_res_wide)), paste0('analysis/output/tables/tab1.tex'))






